 use "$data_out/naive.dta", clear


//reverse-subsample - leave out group and calculate instrument without that group; then apply that instrument to the left out group

*Age < 30
preserve
drop if age < 30
do "$Analysis/Robustness\Physician_Propensity_Reverse"
keep provID year propensity hym hd diag
bys provID year: keep if _n ==1
save "$data_out/propensity_age_old", replace
restore

keep if age < 30
merge m:1 provID year using "$data_out/propensity_age_old", keep(3) nogen
save "$data_out/monotonicity_check_young", replace

*Age >=30
use "$data_out/naive.dta", clear
preserve
drop if age >= 30
do "$Analysis/Robustness/Physician_Propensity_Reverse"
keep provID year propensity hym hd diag
bys provID year: keep if _n ==1
save "$data_out/propensity_age_young", replace
restore

keep if age >= 30
merge m:1 provID year using "$data_out/propensity_age_young", keep(3) nogen
save "$data_out/monotonicity_check_old", replace

*Race - Non_white
use "$data_out/naive.dta", clear
preserve
drop if race_white == 0
do "$Analysis/Robustness/Physician_Propensity_Reverse"
keep provID year propensity hym hd diag
bys provID year: keep if _n ==1
save "$data_out/propensity_white", replace
restore

keep if race_white == 0
merge m:1 provID year using "$data_out/propensity_white", keep(3) nogen
save "$data_out/monotonicity_check_non_white", replace

*Race - White
use "$data_out/naive.dta", clear
preserve
drop if race_white == 1
do "$Analysis/Robustness/Physician_Propensity_Reverse"
keep provID year propensity hym hd diag
bys provID year: keep if _n ==1
save "$data_out/propensity_non_white", replace
restore

keep if race_white == 1
merge m:1 provID year using "$data_out/propensity_non_white", keep(3) nogen
save "$data_out/monotonicity_check_white", replace

*Gender - Male
use "$data_out/naive.dta", clear
preserve
drop if female == 0
do "$Analysis/Robustness/Physician_Propensity_Reverse"
keep provID year propensity hym hd diag
bys provID year: keep if _n ==1
save "$data_out/propensity_female", replace
restore

keep if female == 0
merge m:1 provID year using "$data_out/propensity_female", keep(3) nogen
save "$data_out/monotonicity_check_male", replace

*Gender - Female
use "$data_out/naive.dta", clear
preserve
drop if female == 1
do "$Analysis/Robustness/Physician_Propensity_Reverse"
keep provID year propensity hym hd diag
bys provID year: keep if _n ==1
save "$data_out/propensity_male", replace
restore

keep if female == 1
merge m:1 provID year using "$data_out/propensity_male", keep(3) nogen
save "$data_out/monotonicity_check_female", replace

*Single
use "$data_out/naive.dta", clear
preserve
drop if married == 0
do "$Analysis/Robustness/Physician_Propensity_Reverse"
keep provID year propensity hym hd diag
bys provID year: keep if _n ==1
save "$data_out/propensity_married", replace
restore

keep if married == 0
merge m:1 provID year using "$data_out/propensity_married", keep(3) nogen
save "$data_out/monotonicity_check_single", replace

*Married
use "$data_out/naive.dta", clear
preserve
drop if married == 1
do "$Analysis/Robustness/Physician_Propensity_Reverse"
keep provID year propensity hym hd diag
bys provID year: keep if _n ==1
save "$data_out/propensity_single", replace
restore

keep if married == 1
merge m:1 provID year using "$data_out/propensity_single", keep(3) nogen
save "$data_out/monotonicity_check_married", replace

*Less than College
use "$data_out/naive.dta", clear
preserve
drop if college == 0
do "$Analysis/Robustness/Physician_Propensity_Reverse"
keep provID year propensity hym hd diag
bys provID year: keep if _n ==1
save "$data_out/propensity_college", replace
restore

keep if college == 0
merge m:1 provID year using "$data_out/propensity_college", keep(3) nogen
save "$data_out/monotonicity_check_no_college", replace

*College
use "$data_out/naive.dta", clear
preserve
drop if college == 1
do "$Analysis/Robustness/Physician_Propensity_Reverse"
keep provID year propensity hym hd diag
bys provID year: keep if _n ==1
save "$data_out/propensity_no_college", replace
restore

keep if college == 1
merge m:1 provID year using "$data_out/propensity_no_college", keep(3) nogen
save "$data_out/monotonicity_check_college", replace

*Depression

use "$data_out/naive.dta", clear
merge 1:1 pid_pde encdate using "$data_out/depression_scripts", keep(1 3) gen(depression) keepus(pid_pde encdate)
recode depression 1=0 3=1
save "$data_out/naive_depression", replace

use "$data_out/naive_depression.dta", clear
preserve
drop if depression == 1
do "$Analysis/Robustness/Physician_Propensity_Reverse"
keep provID year propensity hym hd diag
bys provID year: keep if _n ==1
save "$data_out/propensity_no_depression", replace
restore

keep if depression == 1
merge m:1 provID year using "$data_out/propensity_no_depression", keep(3) nogen
save "$data_out/monotonicity_check_depression", replace

use "$data_out/naive_depression.dta", clear
preserve
drop if depression == 0
do "$Analysis/Robustness/Physician_Propensity_Reverse"
keep provID year propensity hym hd diag
bys provID year: keep if _n ==1
save "$data_out/propensity_depression", replace
restore

keep if depression == 0
merge m:1 provID year using "$data_out/propensity_depression", keep(3) nogen
save "$data_out/monotonicity_check_no_depression", replace


*Above Median AFQT_P
use "$data_out/naive.dta", clear
qui summ afqt_p, d
local score `r(p50)'
preserve
drop if afqt_p > `score'
do "$Analysis/Robustness/Physician_Propensity_Reverse"
keep provID year propensity hym hd diag
bys provID year: keep if _n ==1
save "$data_out/propensity_afqt_low", replace
restore

keep if afqt_p > `score'
merge m:1 provID year using "$data_out/propensity_afqt_low", keep(3) nogen
save "$data_out/monotonicity_check_afqt_high", replace


*At or Below Median AFQT_P
use "$data_out/naive.dta", clear

preserve
drop if afqt_p <= `score'
do "$Analysis/Robustness/Physician_Propensity_Reverse"
keep provID year propensity hym hd diag
bys provID year: keep if _n ==1
save "$data_out/propensity_afqt_high", replace
restore

keep if afqt_p <= `score'
merge m:1 provID year using "$data_out/propensity_afqt_high", keep(3) nogen
save "$data_out/monotonicity_check_afqt_low", replace







//Regress and store results


 use "$data_out/provider_propensity.dta", clear

estimates clear


local controls_age "race_white female  junior_enlisted college married longevity afqt_p" 
local controls_afqt "i.age_bin race_white female  junior_enlisted college married longevity" 

local controls_race "i.age_bin  female junior_enlisted married college longevity afqt_p" 

local controls_depression "i.age_bin race_white female  junior_enlisted college married longevity afqt_p" 
local controls_no_depression "i.age_bin race_white female  junior_enlisted college married longevity afqt_p" 
local controls_drugs "i.age_bin race_white female  junior_enlisted college married longevity afqt_p" 
local controls_no_drugs "i.age_bin race_white female  junior_enlisted college married longevity afqt_p" 


local controls_female "i.age_bin race_white junior_enlisted college married longevity afqt_p"
local controls_male "i.age_bin race_white junior_enlisted college married longevity afqt_p"
local controls_single "i.age_bin race_white female junior_enlisted college longevity afqt_p" 
local controls_married "i.age_bin race_white female junior_enlisted college longevity afqt_p" 
local controls_college "i.age_bin race_white female  junior_enlisted married longevity afqt_p" 
local controls_no_college "i.age_bin race_white female junior_enlisted married longevity afqt_p" 






//age < 30
ivreghdfe day_180 `controls_age' (ed_prescription = propensity) if age < 30, absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(st)
use "$data_out/monotonicity_check_young", clear
ivreghdfe day_180 `controls_age' (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(str)
estadd ysumm, replace


label var propensity "Age < 30"

#delimit ;

esttab st*  using "monotonicity.rtf" , 
	replace label se star(* 0.10 ** 0.05 *** 0.01)
	b(4) se(4)
	stats( ,
	fmt(a2 %-9.0fc))
	nonotes
	mtitles("Instrument" "Reverse Sample")
		keep(propensity)
		title(First Stage Results Sub-Sample & Reverse Sub-Sample);

#delimit cr
 

//Age >= 30
use "$data_out/provider_propensity.dta", clear
ivreghdfe day_180 `controls_age' (ed_prescription = propensity) if age >=30, absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(st)
use "$data_out/monotonicity_check_old", clear
ivreghdfe day_180 `controls_age' (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(str)
estadd ysumm, replace

label var propensity "Age > 30"
esttab st*  using "monotonicity.rtf" , ///
	append label se star(* 0.10 ** 0.05 *** 0.01) ///
	b(4) se(4) ///
	stats( ,	fmt(a2 %-9.0fc)) nonotes nomtitles	nonum ///
		keep(propensity) 
		
//Below Median AFQT
use "$data_out/provider_propensity.dta", clear
summ afqt_p, d
local score `r(p50)'
ivreghdfe day_180 `controls_afqt' (ed_prescription = propensity) if afqt_p <= `score', absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(st)
use "$data_out/monotonicity_check_afqt_low", clear
ivreghdfe day_180 `controls_afqt' (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(str)
estadd ysumm, replace


label var propensity "Below Median Qualification Score"
esttab st*  using "monotonicity.rtf" , ///
	append label se star(* 0.10 ** 0.05 *** 0.01) ///
	b(4) se(4) ///
	stats( ,	fmt(a2 %-9.0fc)) nonotes nomtitles	nonum ///
		keep(propensity) 

//Above Median AFQT
use "$data_out/provider_propensity.dta", clear
summ afqt_p, d
local score `r(p50)'
ivreghdfe day_180 `controls_afqt' (ed_prescription = propensity) if afqt_p > `score', absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(st)
use "$data_out/monotonicity_check_afqt_high", clear
ivreghdfe day_180 `controls_afqt' (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(str)
estadd ysumm, replace


label var propensity "Above Median Qualification Score"
esttab st*  using "monotonicity.rtf" , ///
	append label se star(* 0.10 ** 0.05 *** 0.01) ///
	b(4) se(4) ///
	stats( ,	fmt(a2 %-9.0fc)) nonotes nomtitles	nonum ///
		keep(propensity) 

//gender

use "$data_out/provider_propensity.dta", clear
ivreghdfe day_180 `controls_female' (ed_prescription = propensity) if female == 1, absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(st)
use "$data_out/monotonicity_check_female", clear
ivreghdfe day_180 `controls_female' (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(str)
estadd ysumm, replace

label var propensity "Female"
esttab st*  using "monotonicity.rtf" , ///
	append label se star(* 0.10 ** 0.05 *** 0.01) ///
	b(4) se(4) ///
	stats( ,	fmt(a2 %-9.0fc)) nonotes nomtitles	nonum ///
		keep(propensity) 
		

use "$data_out/provider_propensity.dta", clear
ivreghdfe day_180 `controls_female' (ed_prescription = propensity) if female == 0, absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(st)
use "$data_out/monotonicity_check_male", clear
ivreghdfe day_180 `controls_female' (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(str)
estadd ysumm, replace

label var propensity "Male"
esttab st*  using "monotonicity.rtf" , ///
	append label se star(* 0.10 ** 0.05 *** 0.01) ///
	b(4) se(4) ///
	stats( ,	fmt(a2 %-9.0fc)) nonotes nomtitles	nonum ///
		keep(propensity) 
		
// Race

use "$data_out/provider_propensity.dta", clear
ivreghdfe day_180 `controls_race' (ed_prescription = propensity) if race_w == 1, absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(st)
use "$data_out/monotonicity_check_white", clear
ivreghdfe day_180 `controls_race' (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(str)
estadd ysumm, replace

label var propensity "Race - White"
esttab st*  using "monotonicity.rtf" , ///
	append label se star(* 0.10 ** 0.05 *** 0.01) ///
	b(4) se(4) ///
	stats( ,	fmt(a2 %-9.0fc)) nonotes nomtitles	nonum ///
		keep(propensity) 
		

use "$data_out/provider_propensity.dta", clear
ivreghdfe day_180 `controls_race' (ed_prescription = propensity) if race_w == 0, absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(st)
use "$data_out/monotonicity_check_non_white", clear
ivreghdfe day_180 `controls_race' (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(str)
estadd ysumm, replace

label var propensity "Race - Non White"
esttab st*  using "monotonicity.rtf" , ///
	append label se star(* 0.10 ** 0.05 *** 0.01) ///
	b(4) se(4) ///
	stats( ,	fmt(a2 %-9.0fc)) nonotes nomtitles	nonum ///
		keep(propensity) 
		
//Married

use "$data_out/provider_propensity.dta", clear
ivreghdfe day_180 `controls_married' (ed_prescription = propensity) if married == 1, absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(st)
use "$data_out/monotonicity_check_married", clear
ivreghdfe day_180 `controls_married' (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(str)
estadd ysumm, replace

label var propensity "Married"
esttab st*  using "monotonicity.rtf" , ///
	append label se star(* 0.10 ** 0.05 *** 0.01) ///
	b(4) se(4) ///
	stats( ,	fmt(a2 %-9.0fc)) nonotes nomtitles	nonum ///
		keep(propensity) 
		

use "$data_out/provider_propensity.dta", clear
ivreghdfe day_180 `controls_married' (ed_prescription = propensity) if married == 0, absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(st)
use "$data_out/monotonicity_check_single", clear
ivreghdfe day_180 `controls_married' (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(str)
estadd ysumm, replace

label var propensity "Single"
esttab st*  using "monotonicity.rtf" , ///
	append label se star(* 0.10 ** 0.05 *** 0.01) ///
	b(4) se(4) ///
	stats( ,	fmt(a2 %-9.0fc)) nonotes nomtitles	nonum ///
		keep(propensity) 
		
		//College

use "$data_out/provider_propensity.dta", clear
ivreghdfe day_180 `controls_college' (ed_prescription = propensity) if college == 0, absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(st)
use "$data_out/monotonicity_check_no_college", clear
ivreghdfe day_180 `controls_no_college' (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(str)
estadd ysumm, replace

label var propensity "Education < College"
esttab st*  using "monotonicity.rtf" , ///
	append label se star(* 0.10 ** 0.05 *** 0.01) ///
	b(4) se(4) ///
	stats( ,	fmt(a2 %-9.0fc)) nonotes nomtitles	nonum ///
		keep(propensity) 
		

use "$data_out/provider_propensity.dta", clear
ivreghdfe day_180 `controls_college' (ed_prescription = propensity) if college == 1, absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(st)
use "$data_out/monotonicity_check_college", clear
ivreghdfe day_180 `controls_college' (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(str)
estadd ysumm, replace

label var propensity "Education - College Degree"
esttab st*  using "monotonicity.rtf" , ///
	append label se star(* 0.10 ** 0.05 *** 0.01) ///
	b(4) se(4) ///
	stats( ,	fmt(a2 %-9.0fc)) nonotes nomtitles	nonum ///
		keep(propensity) 

		
// Depression		

use "$data_out/provider_propensity.dta", clear
merge 1:1 pid_pde encdate using "$data_out/depression_scripts", keep(1 3) nogen
recode subsample .=0
rename subsample depression

ivreghdfe day_180 `controls_depression' (ed_prescription = propensity) if depression == 0, absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(st)
use "$data_out/monotonicity_check_no_depression", clear
ivreghdfe day_180 `controls_depression' (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(str)
estadd ysumm, replace

label var propensity "Diagnosis - No Depression"
esttab st*  using "monotonicity.rtf" , ///
	append label se star(* 0.10 ** 0.05 *** 0.01) ///
	b(4) se(4) ///
	stats( ,	fmt(a2 %-9.0fc)) nonotes nomtitles	nonum ///
		keep(propensity) 
		

use "$data_out/provider_propensity.dta", clear
merge 1:1 pid_pde encdate using "$data_out/depression_scripts", keep(1 3) nogen
recode subsample .=0
rename subsample depression

ivreghdfe day_180 `controls_depression' (ed_prescription = propensity) if depression == 1, absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(st)
use "$data_out/monotonicity_check_depression", clear
ivreghdfe day_180 `controls_college' (ed_prescription = propensity), absorb(hym hd i.diag mos) cluster(provID) first savefirst savefprefix(str)
estadd ysumm, replace

label var propensity "Diagnosis - Depression"
esttab st*  using "monotonicity.rtf" , ///
	append label se star(* 0.10 ** 0.05 *** 0.01) ///
	b(4) se(4) ///
	stats( ,	fmt(a2 %-9.0fc)) nonotes nomtitles	nonum ///
		keep(propensity) 

				
		

